Force distributions near the jamming and glass transitions 
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We calculate the distribution of interparticle normal forces P{F) near the glass and jamming 
transitions in model supercooled liquids and foams, respectively. P(F) develops a peak that ap- 
pears near the glass or jamming transitions, whose height increases with decreasing temperature, 
decreasing shear stress and increasing packing density. A similar shape of P(F) was observed in 
experiments on static granular packings. We propose that the appearance of this peak signals the 
development of a yield stress. The sensitivity of the peak to temperature, shear stress and density 
lends credence to the recently proposed generalized jamming phase diagram. 
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Granular materials can flow when shaken, but jam 
when the shaking intensity is lowered (!]]. Similarly, 
foams and emulsions can flow when sheared, but jam 
when shear stress is lowered || . These systems are ather- 
mal because thermal energy is insufficient to change the 
packing of grains, bubbles, or droplets. When the exter- 
nal driving force is too small to cause particle rearrange- 
ments, these materials become amorphous solids and de- 
velop a yield stress. A supercooled liquid, on the other 
hand, is a thermal system that turns, as temperature is 
lowered, into a glass — an amorphous solid with a yield 
stress H . Despite significant differences between driven, 
athermal systems and quiescent, thermal ones, it has 
been suggested that the process of jamming — developing 
a yield stress in an amorphous state — may lead to com- 
mon behavior, and that these systems can be unified by 
a jamming phase diagram 0]. This implies that there 
should be similarities in these different systems as they 
approach jamming or glass transitions. We test this spec- 
ulation by measuring the distribution P(F) of interpar- 
ticle normal forces F, in model supercooled liquids and 
foams. We find that for glasses, P(F) is quantitatively 
similar to experimental results on granular materials || . 

When granular materials jam, the distribution of 
stresses is known to be inhomogeneous ^|J7) . As proposed 
in Ref. M, we quantify this effect by measuring P(F). 
Our aim is to determine which feature in P(F) is asso- 
ciated with development of a yield stress. Experiments 
HU and simulations |],[l(J on static granular packings 
find that P(F) has a plateau or small peak at small F 
and decays exponentially at large F. We argue that the 
development of a peak is the signature of jamming. 

For supercooled liquids, equilibrium statistical me- 
chanics gives insight into the shape of P(F). Since forces 
depend only on particle separations, P(F)dF = G(r)dr, 
where G(r)dr is the probability of finding a particle be- 
tween r and r + dr given a particle at the origin. Thus, 
G(r) — p/(N— l)SDr D ~ 1 g(r), where N is the number of 



particles, p is the number density, g(r) is the pair distri- 
bution function, and Sdt 0-1 is the surface area of a D- 
dimensional sphere of radius r. Although it is well known 
that g(r) does not change significantly as the tempera- 
ture is varied through the glass transition T g , we show 
below that P(F) is quite sensitive and actually develops 
a peak near T g . Physically, forces (or stresses JIT| , ^4| ) 
are crucial for understanding the slowing down of stress 
relaxation near the glass transition, or the development 
of a yield stress. It is therefore not surprising that P(F), 
which is a particular weighting of g(r), is much more 
sensitive to the glass transition than g(r) itself. 

In a jammed system like a granular material, an an- 
alytic expression for g(r) is not known and P(F) must 
be measured directly. However, in an equilibrium system 
at temperature T, the large-force behavior of P(F) can 
be obtained from the small-separation behavior of g(r): 
g(r) = y(r) exp[— V(r)/kbT], where V(r) is the pair po- 
tential and y(r) depends relatively weakly on r at small 
r |Ol. This leads to: 



P[F(r)} ~ y(r)r D - 1 %exp[-V(r)/k b T\ 
dh 



(1) 



From our simulations, we compute the force distribu- 
tions in systems that are out of equilibrium, such as 
glasses and sheared foams, as well as systems in equi- 
librium, such as supercooled liquids. We find that P(F) 
for supercooled liquids (with sufficiently strong repulsive 
potentials) decays approximately exponentially at large 
forces, as predicted by Eq. |l|. Because this is true at all 
temperatures, even those far above T g , the exponential 
tail is not necessarily a signature of an amorphous solid. 

We perform constant-temperature molecular dynamics 
simulations on binary mixtures in 2D, using the Gaussian 
constraint thermostat and leapfrog Verlet algorithm [jl3| . 
The masses m of the particles are the same, but the ra- 
tio of particle diameters, <7ij<J\ = 1.4, ensures that the 
system does not crystallize Jul. We confine N — 1024 
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particles (512 of each variety) to a square box and use 
periodic boundary conditions. For each simulation, we 
choose one of the following interparticle pair potentials: 



V^{r)=e{a ab /r) 
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(2) 



V^(r)=4e {a ab /r)^-{a ab /rf 



V a L b ,m2 (r) ee 4e [K fc /r) 12 - (a ab /rf] + e; r/a ab < 2 
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where daf, = (a a + <J b )/2 for a,b — 1,2. The potentials 
LJR12 and LJR2A are zero above the specified cutoffs. 
(The potentials SC and LJ are truncated at large r, 
r/cah = 4.5.) Below, we measure time, force, temper- 
ature, and density in units of o~i{m/ 'e) 1 / 2 , e/ci, e/fcb, 
and erf 2 respectively. The simulations on purely repul- 
sive potentials SC, LJRY2, and LJR2A simulations were 
carried out at constant density p = 0.747; the simula- 
tions on LJ, which includes an attraction, were carried 
out at zero average pressure. 

The hallmark of the glass transition is the extreme 
slowing down of the dynamics as temperature is lowered 
toward the glass transition. The pair potentials in Eq. ^ 
all give rise to glass transitions as temperature is lowered 
fuj . We determine T g by measuring the self-part of the 
intermediate scattering function F2(k p ,t) for the large 
particles at a wavevector k p corresponding to the first 
peak of the static structure factor Jl4],[l5| • For high tem- 
peratures, the liquid equilibrates quickly and F2{k p ,t) 
decays exponentially to zero. The relaxation time r r is 
defined as the time at which F2{k p ,t) decays to 1/e; this 
is a measure of the a-relaxation time |p^ , |l5| . Since r r 
increases so rapidly near the glass transition, simulations 
can only reach equilibrium for temperatures T > T gi 
where T g is determined by when r r exceeds a prede- 
termined, large value, which we take to be r r > 1000. 
For our parameter choices, the glass transition tempera- 
tures are T g sc = 0.38, T s iJm2 = 1.1, T g LJim = 3.0 and 
T 9 iJ = 0.17. 

For T > T g we measure P(F) for all interparticle force 
pairs from at least 250 configurations after equilibrat- 
ing each configuration for 10 — 100r r . The top frame 
of Fig. |l| shows P(F) plotted versus F/T for LJR12 for 
seven temperatures above T g with temperature decreas- 
ing from top to bottom. At high temperatures, we see in 
the top frame of Fig. [I] that P(F) increases with decreas- 
ing F over the entire range of F. However, as tempera- 
ture is lowered towards T g , a plateau in P(F) forms at 
forces below the average (F) . By shifting each curve ver- 
tically, we have obtained collapse of the high-force data 
for all of these equilibrium systems. From Eq. ^, we 
expect that the large-force tail should scale asymptoti- 
cally as exp(— BF 12 / 13 /T), where B is a constant and 
the power 12/13 derives from the 1/r 12 repulsion. Thus, 



for particles with harder cores (steeper repulsions), the 
tail becomes closer to an exponential in F, as seen in 
experiments on granular materials [E| . This explanation 
for the exponential tail is different from that of the q- 
model J?J and its generalizations [|l6| based on stochas- 
tic force propagation. Previous LJ simulations along the 
liquid-vapor coexistence line |l7j showed that the Carte- 
sian components of the force also have an exponential 
distribution. Our results are related to theirs: for high 
forces, the total force on a particle, which is the vector 
sum of the normal forces, will be dominated by the largest 
normal force. This is why the distribution of Cartesian 
components is also exponential. 




FIG. 1. Top: P(F) for all interparticle force pairs versus 
F/T for LJR12 obtained for seven temperatures above T g with 
T decreasing from top to bottom. Bottom: P(F/(F}) versus 
F/(F) for LJR12 for two temperatures above and three below 
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We also study P(F) out of equilibrium by performing 



thermal quenches from Tj > T g to T i 



<T g . 



The results 



discussed below are relatively insensitive to changes in Tj 
or quench rate. In the bottom frame of Fig. ^, we show 
the long-time behavior of P(F) for LJR12 following a 
quench below T g w 1.1 to T f = 0.6, 0.3, and 0.1. For 
comparison, two equilibrium distributions at T — 1.5 and 
3.0 are also shown. We scaled the abscissa by (F) (which 
increases with T). There are two significant features in 
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P(F) for glasses. First, the slope of the exponential tail 
increases as T is lowered. The temperature correspond- 
ing to the tail T ta ii, however, is not the final temperature 
Tf, but rather satisfies Tf < T ta y < T g . Thus, a frac- 
tion of the large thermal forces cannot relax in the glassy 
state. The second significant feature of P(F) for glasses is 
the formation of a peak near (F), as shown in the bottom 
frame of Fig. [j]. Thus, in contrast to g(r), there is a sig- 
nificant change in P(F) below T g . The behavior of P(F) 
when quenched below T g is qualitatively the same for all 
potentials and densities studied, showing that the peak 
signals the glass transition in a system with attractive 
interactions and no applied pressure as well as systems 
with purely repulsive interactions under pressure. 
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FIG. 2. P(F/{F)) versus F/(F) for both LJR12 and LJR24 
after a quench to Tf = 0.8. Data from experiments on static 
granular packings from [^| are also shown. Note that all three 
sets of data have a plateau at small F and decay exponentially 
at large F. 



The potentials LJR12 and LJR24 in Eq. || are most 
similar to granular materials since they produce purely 
repulsive forces that vanish at small separation. We com- 
pare P(F/(F)) in the glassy state for LJFLYl and LJR24 
to P(F/(F)) for static granular packings in Fig. ||. For 
both LJR12 and LJR24, we have quenched to Tf = 
0.8 < T g . Remarkably, the force distributions, when 
scaled by the average force (F), are nearly identical 
for LJR12, LJR24, and experiments on static granular 
packings j^] over the entire range of forces. This implies 
that for sufficiently hard repulsive potentials, the shape 
of the distribution is not sensitive to the shape of the 
potential. In the limit of hard spheres, where the power 
of the repulsive term in the potential diverges, we ex- 
pect similar behavior. In systems with softer potentials, 
such as Hertzian or harmonic repulsive springs, we also 
find the same shape of P(F) as in Fig. ^ at very low 
temperatures near the close- packing density Ej . These 
results suggest that the slight peak or plateau at small 
forces and exponential tail at large forces are a generic 
feature of P{F) in athermal, experimental systems near 
the onset of jamming. 
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3. Top: P{F/{F)) versus F/(F) for foams with 
for several (j) near random close-packing. Bottom: 

xy 



P(F/{F)) for foams with (p — 0.9 > 4>o and a xv lowered to- 
wards a 



ti- 



ls a peak or plateau in P(F) also observed in other 
jammed systems? To answer this, we have studied model 
two-dimensional foams |18jl^], where bubbles are treated 
as circles that can overlap and interact via two types 
of pairwise interactions. The first is a harmonic repul- 
sion that is nonzero when the distance between centers 
of two bubbles is less than the sum of their radii. The 
other is a simple dynamical friction proportional to the 
relative velocities of two neighboring bubbles. In foam, 
thermal motion of bubbles is negligible. We simulate a 
400-bubble system at constant area with periodic bound- 
ary conditions in the ^-direction and fixed boundaries in 
the y-direction. Bubble radii Ri are chosen from a flat 
distribution with 0.2 < R l /(R) < 1.8. 

At packing fractions <p above random close packing 
((fro ~ 0.84), quiescent foam is an amorphous solid with 
a yield stress a y . However, when shear stress a xy > a y 
is applied, the foam flows. There are therefore two ways 
to approach the amorphous solid. We can either increase 
<p towards cpo at a xy = (route 1), or we can decrease 
a xy towards a y at fixed <p > (po (route 2). In Fig. we 
show P(F/(F)) (only including harmonic elastic forces) 
along these two routes. The distributions along route 1 
in the top frame were measured after quenching 50 con- 
figurations from <pi <C <po to <j> by increasing each particle 
radius. When <p < <po, P(F/(F)) increases monotoni- 
cally as F/(F) decreases. As <p increases above <po, a 
local maximum forms near (F) . A similar trend is found 
along route 2. To obtain these distributions, we aver- 
aged over at least 500 configurations with each brought to 
steady state for a strain of « 10. In all cases shown, a xy 
exceeds cr y , so the systems are flowing. We find that at 
large a xy , P(F/ (F)) is nearly constant at small F. When 
a xy is lowered towards a y w 0.10, a peak in P(F/(F)) 
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forms near F/(F) ~ 1. Similar behavior is observed in 
P(F/(F)) as a function of in experiments on sheared 
deformable disks |2(J and as a function of confining stress 
in simulations of deformable spheres JTc| ]. 

In this Letter, we have shown a connection between 
development of a yield stress, either by a glass transition 
or conventional jamming transition, and the appearance 
of a peak in P(F). We have established that four dif- 
ferent model supercooled liquids develop first a plateau 
and then a peak in P(F) as temperature is lowered below 
the glass transition. We have also found (but not shown 
here) that the LJR12 liquid displays identical results for 
P(F) as shear stress is lowered from the flowing state or 
as density is raised from the liquid state at fixed temper- 
ature |2lj . We also find that the appearance of a peak in 
P(F) coincides with onset of crystallization in monodis- 
perse liquids, further supporting the connection between 
the peak and the yield stress plfl . The athermal foam 
likewise develops a peak in P(F) as it approaches jam- 
ming along two different routes. Static granular packings 
exhibit a plateau or small peak in P(F) as well. Thus, a 
peak in P{F) appears as a wide variety of systems jam 
along each of the axes of the jamming phase diagram [Q . 
This suggests that jamming leads to common behavior 
and that the glass transition may resemble more conven- 
tional jamming transitions. 

This still leaves the question of why formation of a peak 
or plateau in P(F) appears to signal the development of a 
yield stress. The presence of the peak or plateau implies 
that there are a large number of forces near the aver- 
age value. This is consistent with the existence of force 
chains, since each particle within a force chain must have 
roughly balanced forces on either side. We speculate that 
systems jam when there are enough particles in a force 
chain network to support the stress over the time scale 
of the measurement. Forces at the peak of P(F) are 
among the slowest to relax: these forces correspond to 
separations near the first peak of g(r), which stem from 
wavevectors near the first peak in the static structure 
factor, which are among the most slowly relaxing modes 
p2| . This implies that force chains observed in granular 
packings may also be important to the glass transition. 
The fact that force chains do not couple strongly to den- 
sity fluctuations may explain why they have not been 
observed directly. However, large kinetic heterogeneities 
that appear near T g ||| may be linked to the formation 
of force chains. This interpretation suggests that force 
chains may provide the key to the elusive order parame- 
ter for the glass transition. 
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